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We report on electronic transport calculations for self- assembled mono-layers (SAM) of 1,4- 
phenylene diisocyanide on Au(lll) contacts. Experimentally one observes more structure (i.e peaks) 
within the measured conductance curve for this molecule with two cyanide end-groups, compared 
to measurements with molecules having thiol end-groups. The calculations are performed on the 
semi-empiric extended Hiickel level using elastic scattering quantum chemistry (ESQC) and we in- 
vestigate three possible explanations for the experimental findings. Comparing the experimental 
and theoretical data, we are able to rule out all but one of the scenarios. The observed additional 
peaks are found to be only reproduced by a mono-layer with additional molecules perturbing the 
periodicity. It is conjectured that the weaker coupling to Au of cyanide end-groups compared to 
thiol end-groups might be responsible for such perturbations. 

PACS numbers: 85.65.+h, 73.23.-b, 72.10.-d 



I. INTRODUCTION 

Within the last decade an increasing interest in molec- 
ular electronics has developed, with the expectation of 
realising molecular diodes and transistors. This is based 
on the progress in manipulation techniques, which now 
allow the controlled attachment of atomic scale struc- 
tures like molecules to mesoscopic leads. With these new 
devices one is able to determine the conductance prop- 
erties of molecular structures. Explaining and predicting 
the electronic behaviour of such devices is an essential 
step towards their design and use as nano-scale electronic 
circuits. 

To this end a number of theoretical studies have been 
performed with the aim of reproducing measured IV 
characteristics. These studies differ in the way they take 
the electronic levels of the molecules, their modification 
by the coupling to the leads and the change of electro- 
static potential due to bias into account. Semi-empiric 
methodsi*^ have been used, as well as first principles 
techniques 4*^^ the latter being restricted to systems of 
moderate size. 

The wide range of experimentally observed behaviour 
(see section ITT|) suggests, that not only the structure of 
the molecule, but also the details of the device fabrica- 
tion affect the conduction properties of molecular devices. 
The crucial step is the deposition of molecules onto the 
surface of the lead. As this is achieved by self-assembly 
the amount of adsorbed molecules and their individual 
positions can not be exactly controlled and therefore re- 
mains unknown. A satisfactory understanding of the in- 
terplay between geometrical alignment of the molecules 
and measured conductance properties has thus not yet 
been achieved (for a recent review see e.g. Nitzan et al£). 

In this paper we study the way in which changing the 
geometrical alignment of the mono-layer has an influence 
on the conduction properties of a molecular device. In 
so doing we can rule out a number of explanations which 



have previously been considered^ to explain the occur- 
rence of additional structure in the conductance-voltage 
(CV) characteristics. 

The outline is the following: first we summarise some 
of the recent experimental findings. Then the method we 
use (based on ESQC 10 ) for calculating the conductance 
properties of molecular devices, is discussed. Calcula- 
tions for the conductance of a SAM, being build of 1,4- 
phenylene diisocyanide (PDI) and sandwiched between 
gold leads are then presented. The results for three qual- 
itatively different geometrical constitutions of the mono- 
layer are compared with experimental data. By this we 
can rule out all but one and conclude that the only ge- 
ometrical alignment, which gives rise to several peaks in 
the conductance curve, is a mono-layer with additional 
molecules perturbing the periodicity. 



II. EXPERIMENTAL OVERVIEW 

The devices built to study conductance properties of 
molecular structures differ not only in amount and chem- 
ical structure of the molecules in use but also in the way 
these are attached to metallic or semi-conducting leads. 
Single or few molecules are accessible in mechanically 
controllable break junctions (MCB) and with the scan- 
ning tunneling microscope (STM). Many molecules are 
involved in sandwiched self-assembled mono-layer (SAM) 
experiments. The observed properties depend on the ex- 
act geometry of the device. The conductance differs in 
orders of magnitude and the qualitative voltage depen- 
dency of the current ranges from simple ohmic behaviour 
to negative differential resistance (NDR). 11 

In the past Reed et have measured the electri- 
cal conductance of a self-assembled molecular mono-layer 
bridging an MCB at room temperature. Molecules of 1,4- 
benzene dithiol (i.e. having two thiol end-groups, which 
are known to couple strongly to Au-atoms) were used 
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and the CV characteristic was found to be symmetric 
with one peak in the voltage range of — 2V. They mea- 
sured a current of the order of 50nA at a bias voltage of 
2V, which they claim is produced by transport through 
one single active molecule. Reichert et al^ also used an 
MCB with molecules having two thiol end groups, but 
being considerably longer. The measured current ampli- 
tude was about 500nA at IV, i.e. although the molecule 
was more than twice as long, the current was ten times 
larger. 

With a different setup, where a SAM is sandwiched be- 
tween two metallic leads, Chen et alii have found nega- 
tive differential resistance (NDR) , namely one peak at 2V 
in the IV curve. The molecule under investigation had 
one thiol end-group only and was attached to Au-leads 
at both ends. The measurements were taken at room 
temperature and the measured current maximum was of 
the order of InA. 

Only recently, sandwiched SAM devices at 4.2K were 
studied(2ii4 where a benzene ring with two cyanide in- 
stead of thiol end-groups was used. The measurements 
exhibited currents of the order of 50 — 400nA. The CV 
characteristic for this molecule revealed more structure, 
in form of three to five peaks within a voltage range of 
IV. Such a behaviour was not observed with previous 
devices containing other molecules. 



III. THEORETICAL FORMALISM 

In the literature there has been presented quite a num- 
ber of techniques to calculate non-equilibrium electronic 
transport through molecular systems attached to meso- 
scopic leads. Usually the Landauer formalism is applied, 
which describes current as elastic electron transmission 
and therefore requires the transmission function T(E). 
To this end one needs a framework that allows for a 
description of the molecular device on the atomic level. 
This involves not only the molecules themselves, but also 
the surface and bulk region of the leads. Quantum chem- 
istry provides such a framework, and one can choose the 
level of theory according to the size of the system under 
consideration and the computational effort one is willing 
to spend. 

Using a quantum chemistry method, the transmission 
function can be obtained from an effective one-particle 
Hamiltonian, which is an appropriate description for 
strong coupling of the molecules to the leads (as in the 
case of covalent binding) . The methods differ in the gen- 
eration of the one-particle Hamiltonian, which might be 
based on semi-empirical groundsill or on first-principles 
and self-consistent techniquesA^Sii 

A different approach^ taking many-particle effects ex- 
plicitly into account, uses a master equation with transi- 
tion rates calculated perturbatively using the golden rule. 
This approach is appropriate for weak coupling. 

We use the Landauer formalism, as the molecules are 
assumed to be chemically bonded to the gold contacts 



(i.e. strong coupling), together with the semi-empiric ex- 
tended Huckel quantum-chemistry method ESQC. 10 The 
molecular structure is optimised^ beforehand. This ap- 
proach, though limited as compared to more sophisti- 
cated quantum chemistry methods, is yet justified, be- 
cause we want to gain a qualitative understanding of a 
many molecule experiment, which can not be described 
by first-principle techniques, as the number of atoms in- 
volved is beyond the practical limitations of to-date com- 
puter resources. 



A. Landauer formalism 



According to the Landauer formula, current along a de- 
fect region is the result of electron transmission from the 
source to the drain lead. This is described by the trans- 
mission function T(E) and net electron transfer happens 
at all those energies, where there are more states occu- 
pied in the source lead than in the drain. The difference 
in occupation is a result of the applied bias V, which 
raises the chemical potential fix of the source lead above 
the one of the drain lead /Lt2 = (J>i — eV and thus changes 
the distribution function fi(E) = f(E — fii). 



— 2e f°° 

I=—J T(E)(f(E f(E-to))dE. (1) 



The Landauer formula is valid under the condition that 
transport is coherent across the molecule, which is plau- 
sible as the typical mean-free path of electrons within 
metals is of the order of 500nm, while the molecular gap 
between source and drain lead is only about I — 5nm in 
length. Furthermore it is assumed that the coupling of 
the leads to macroscopic reservoirs is reflection-less, and 
that the tunneling electrons equilibrate only deep inside 
the leads. This ensures that the distribution function for 
incoming electrons in a lead can be taken to be spatially 
constant, even close to the molecular region. 

The system is formally partitioned into three regions 
<G {0, 1,2}, two of them (Si, 2) containing the semi- 
infinite leads, the third one (So) being the finite region 
containing all molecules as well as a few surface layers of 
each lead (see Fig. ^| . We use periodic boundary condi- 
tions in the directions perpendicular to the surface nor- 
mal. 

By a tight binding approximation, the infinite- 
dimensional Hamiltonian of the entire system can be 
composed of quantum chemical one-particle block Hamil- 
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channels are connected to each other via the molecular 
region. This is described by the scattering matrix and 
the transmission function is finally obtained by summing 
up the contribution from each channel. 

The calculation can be performed either using Green's 
function techniques^ or equivalentlyi^ using ESQGi2iiS ; 
which is a scattering-matrix approach. We present the 
details of the calculation in the second scheme, as individ- 
ual contributions from each channel to the transmission 
function can then be easily studied. 



B. Bulk propagator 




FIG. 1: (Colour online) Partitioning of the system into three 
parts: the two semi-infinite leads Ei^ (surrounded by boxes) 
and the molecular region Eo between them. 
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The first summation describes the isolated molecular re- 
gion, by an on-site energy e and a hopping term. The 
indices i and j run over the orbital basis set within that 
region. The next two summations describe the isolated 
leads, labelled by d. Layer by layer, starting with the 
surface layer Iq, the first term accounts for intra-layer 
interactions, while the second one describes the interac- 
tion between layers. (The size of each layer is chosen 
such that only adjacent layers have non-zero interaction. 
It therefore depends on the details of the tight-binding 
approximation.) Finally the last term describes the cou- 
pling between the molecular region and each lead. Note 
that only the first layer l contributes to that term and 
that there is no interaction between different leads. These 
are only formal restrictions, as parts of each lead can be 
included into the molecular region. 

The determination of the transmission function in- 
volves two steps. First the conduction properties of the 
isolated leads have to be calculated. Thereby each lead 
will be decomposed into conducting and non-conducting 
in- and out-going channels. These correspond to prop- 
agating and evanescent solutions moving in one of two 
possible directions respectively. In a second step, the 



First we will restrict our attention to the semi-infinite 
lead Hamiltonians, which do not have to be identical. 
The Hamiltonians of equation @ for one lead Hf v are 
layer independent, if one assumes periodicity, i.e. Hfi = 

H ?oh and H W+i = H ti a +i- Usin S Bloch's theorem one 
can reduce the infinite dimensional system of equations 
to an JV x iV-matrix equation (N being the number of 
orbital basis functions in one layer) 

(M(E) + h(E)e tkA + h' ! (E)e- ikA ) 7 (Jfe, E) = 0, (3) 



with M(E) := H h i 
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ESi j +i, and Sw is the overlap matrix between orbitals 
in layer I and layer I' for cases when one does not deal 
with an orthonormal basis set (otherwise Su> = Id ■ <5;;')- 
With A we denote the lattice spacing. Defining A 



one can easily see that equation © is an N x N 
quadratic eigenvalue equation. It can be transformed into 
a 2N x 2N linear eigenvalue equation: 



P(E) 
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(where we have dropped the energy dependency for no- 
tational ease) . This layer-to-layer propagator P reduces 
the problem of finding solutions for the entire isolated 
lead Hamiltonian to specifying the wave function coeffi- 
cients at two adjacent layers 7; and 7/+1 only. This is due 
to the fact that given these coefficients any other pair of 
coefficients 7;+2fe, 7/+2fc+i can now be obtained via: 



ll+2k 
7/+2fc+l 
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(6) 



All possible solutions at energy E can be decomposed 
into independent channels, by solving for the eigenvalues 
of equation These eigenvalues come in pairs such 

that for each eigenvalue A>, there exists a corresponding 
eigenvalue A< satisfying the relation A> = 1/A<, as can 
be seen by transposing equation (J3J. Eigenvalues with 
|A| ^ 1, i.e. complex k, belong to exponentially diverging 
solutions (see equation (0). These are of course non 
physical, as long as the lead is infinite. In semi-infinite 
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leads however (which we are dealing with) , exponentially 
decaying coefficients at the boundary will contribute to 
the surface wave function and must not be neglected. 



C. Current operator 

The contribution from a single channel to the net cur- 
rent can not directly be seen from equation It de- 
pends on the current density associated with a solution to 
the Schrodinger equation ihdtS^f — H7 and is obtained 
via the continuity equation. The probability amplitude 
7I 2 for a stationary solution is constant in time 

l^S-y = ~ (7^7 - 7^7)) = 0, (7) 

because H and S are hermitian. For the probability am- 
plitude at all layers between l\ and I2 one therefore has 
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with the projectors (^7) := 7 ; . This gives rise to the 
definition of the current operator Wi for layer I as 



-h 
h) 



(9) 



Now let both <p and & be solutions at fixed energy E 
with the eigenvalues Ai and A2 respectively. Because the 
expectation value for Wi is layer independent (equation 
El one has: 



(#\Wi\(p) =(§\Wi +1 \ip) 

=X* 1 \ 2 ^\Wi\ip). 



(10) 



This equation describes the connection between the cur- 
rent properties of a solution (p and its eigenvalue A. We 
summarise the results of a detailed analysis of this equa- 
tion, which is given in appendix 1X1 Each channel \ipi) 
can be assigned a current value Vi, defined as 



vt := lm(ipi\W\ifii), 



(11) 



where we have used the layer independence of Wi in sim- 
ply writing W. 



Channels with eigenvalue modulus |A| 7^ 1, i.e. evanes- 
cent waves have zero current value. They therefore do 
not contribute to the current. (Yet they are important 
at the surface, as already mentioned above.) 

Only channels with an eigenvalue of modulus one 
(|A| = 1) contribute to the current. The sign of m deter- 
mines the direction of charge transport. 

Therefore the eigenvalues can be sorted into incoming 
and outgoing solutions (A> and A< respectively), accord- 
ing to the following scheme: if |A| < 1 it represents an 
incoming evanescent solution, if |A| > 1 then it is outgo- 
ing evanescent. Only if |A| = 1 it belongs to a propagat- 
ing solution, which is incoming for v > and outgoing 
otherwise. 

We now define A> and A< as the two N x N diagonal 
matrices composed of all incoming and outgoing eigen- 
The 2N x 2iV-matrix U, which 



values := diag(A^) 
diagonalises P: 



A> 
A< 



has the following quadratic block form: 



U> U < 
J7>A> [/<A< 



(12) 



(13) 



After this transformation into the diagonal basis of the 
propagator, we can easily obtain all physically relevant 
solutions of the infinite lead by specifying the amplitudes 
of all propagating waves at one lattice site. 



D. Scattering matrix 

Up to now, we have considered the isolated leads only. 
These are now assumed to be each coupled to the molec- 
ular defect region and thereby indirectly coupled to one 
another. We are interested in stationary solutions, which 
consist of an incoming propagating wave in one lead, be- 
ing scattered among all the accessible outgoing channels 
(propagating and evanescent ones). This information is 
contained in the scattering matrix S 



six S12 
S21 s 2 2 



(14) 



which determines the wave amplitudes of all outgoing 
waves B, C given the incoming ones A, T>. In appendix IbI 
it is shown, that S is of the form 
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-Mr 1 ■ M n 



(15) 



It is important to notice, that the scattering matrix 
is always quadratic, because in each lead there are the 
same amount of incoming and outgoing channels. This is 
opposed to the transfer matrix Tt ran sf, which determines 
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FIG. 2: Incoming and outgoing wave amplitudes. 



the amplitudes of in- and outgoing waves C, T> in one lead 
given the in- and outgoing waves A, B of a second lead: 



transf 



(16) 



This matrix is quadratic only if both leads have the same 
number of channels. It then is of the form 20 



-^transf 



F C?t 
G Ft 



and the relation to the scattering matrix is2°« 

S = 



(17) 



(18) 



Methods calculating the scattering matrix via the trans- 
fer matrisi^ fail, if two types of leads are used, because 
F is then no longer quadratic and can not be inverted. 
Therefore one commonly takes source and drain lead to 
be identically constituted. But even in such cases, the 
method becomes numerically unstable, with increasing 
distance between the molecular region and one lead. This 
is because the matrix elements of F and G (in equation 
HI 7( 11 diverge exponentially, with increasing lead separa- 
tion. Taking the inverse of F is therefore a numerically 
critical procedure. Both these problems are avoided by 
the direct calculation of the scattering matrix, which we 
present in appendix [BJ This calculation is well defined 
without any restrictions to the number of leads and their 
composition. That means that it is not necessary to re- 
strict to identical leads. Furthermore it allows a numeri- 
cally stable determination of the scattering matrix, even 
for weak coupling. 



E. Transmission function 

The transmission function is the sum over the contri- 
butions from each combination of incoming channel in 
the source lead to outgoing channel in the drain lead: 
T(E) = J^i j Tj<-i- The relation between scattering ma- 
trix S and these transmission function elements is 



T, 



(S2i)j^| 2 



(19) 



where s 2 i is the lower left block of S as defined in equa- 
tion (|14fl . The weighting with velocity factors comes 



about because the scattering matrix s does not relate 
current densities, but wave amplitudes. The current den- 
sities are obtained from these wave amplitudes by multi- 
plication with the corresponding velocity factor Vj. The 
factor Vi in the denominator normalises the transmis- 
sion function to be exactly one for perfect transmission. 
If Vi = then Tj^. L = 0, because incoming evanescent 
waves have zero amplitude at the surface. 

The total current is made up of the contribution from 
each fc-stateii: 



I = 



^^ViT^iEk) (f(E k - Ml ) - f(E k - M2 )) 
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T(E)Af(E)dE. 



Here the summation over k has been transformed into an 
integral over E with a factor of 2 accounting for spin a. 



IV. CALCULATIONS FOR PDI 

Low temperature experiments with PDI SAM's sand- 
wiched between two metallic leads show several peaks in 
the CV-diagrami*^ The typical voltage differences of 
these peaks are in the range of AU ks 0.2V (i.e. there are 
about 5 peaks within U = 0V and U — IV). The com- 
monly adopted explanation for the occurrence of such 
peaks is the following. Each molecular orbital that en- 
ters the energy window, which is opened by the applied 
voltage, enables resonant tunneling. This increases the 
conductance and therefore results in a peak within the 
CV-diagram. 

Typically, the energy gap between molecular orbitals 
is in the range of AE w leV. In other words, for ap- 
plied voltages up to U = IV there should be only a sin- 
gle accessible orbital per molecule, giving rise to only 
a single peak in the CV-diagram. Therefore the fol- 
lowing question arises: are there geometrical alignments 
of the molecules such that the additional peaks in the 
CV diagram can also be explained by resonant tunneling 
through molecular orbitals? 



A. Influence of changes in the molecular alignment 
to the transmission spectrum 

During the device fabrication, the step under least ex- 
perimental control is the adsorption of the molecules onto 
the leads. Therefore the exact geometrical alignment of 
the molecular SAM and, at least in the sandwich geom- 
etry, also the atomic shape of the top metallic lead, is 
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not exactly known. One therefore has to expect not only 
one specific but rather quite a variety of molecular align- 
ments to be produced. As one is interested in the con- 
duction properties of the resulting device, it is important 
to understand the influence of each type of geometrical 
alignment to the transmission function. 

To this end, we have investigated three such possible 
alignments, which will be discussed separately. 

1. By adding metal- atoms on top of the molecular 
mono-layer in a random way, there might occur 
metallic clusters on top of the mono-layer. These 
affect the electronic configuration of the molecules 
individually, and might therefore have an influence 
to the transmission function. 

2. In a SAM experiment, there is not just one 
molecule, but rather a few hundred molecules in- 
volved. If the contribution to the transmission 
function was different for each molecule, then T(E) 
would change qualitatively, with a change in size 
of the mono-layer. We therefore analyse how the 
transmission functions depends on the number of 
molecules involved. 

3. If the molecular mono-layer is not strictly periodic, 
then there will be defects. For example can the dis- 
tance between two molecules be reduced, such that 
inter-molecular bonds can be build. Each of these 
defects will have a specific electronic structure and 
will therefore influence the transmission function. 



1. Influence of metallic clusters 

In the sandwich geometry, first the bottom metallic 
lead is created. Then the molecular mono-layer is ad- 
sorbed on top of it by self-assembly. Finally the top 
metallic lead is build upon the molecular mono-layer. 
The exact shape of neither metallic surface is known and 
may be anything but flat and regular. 

It is likely that the surface atoms of the top metallic 
lead build up clusters on top of the molecular layer (as 
for example in Fig. 0b) . Which influence do they have on 
the electronic configuration of the molecule they are in 
contact with? And do the clusters act as small molecules 
with new electronic levels? 

The influence of an Au cluster on the molecular elec- 
tronic structure is twofold. First it introduces new elec- 
tronic levels, and second the existing molecular electronic 
levels will be shifted, by an amount which depends on the 
strength of the coupling between cluster and molecule. 

The latter effect will be observed as a shifted peak in 
the transmission function, only if the coupling between 
cluster and molecule is different to the coupling between 
top electrode and molecule. For clusters similar to the 
one shown in Fig. [31(b), this is however not the case. The 
energetic peak positions are identical, as can be seen in 
Fig. 0(c). 



Furthermore, there are no additional peaks, which one 
might have expected because of the additional electronic 
levels of the cluster. The explanation of their absence is 
the following: an electronic level gives rise to a peak in 
the transmission function only, if the corresponding or- 
bital wave function overlaps with both the top and bot- 
tom electrode. The overlap with the electrode the clus- 
ter is attached to (say top electrode) is of course large. 
The overlap with the bottom electrode consists of two 
parts. The direct overlap and the indirect overlap via 
the molecule. The direct overlap is negligible due to the 
large spatial separation. The indirect overlap depends 
on the molecular orbital wave function. If the energy of 
the cluster level does not coincide with a molecular en- 
ergy level, then there is no indirect overlap. Only if two 
levels coincide, the indirect coupling is large, but in that 
case, there already exists a transmission peak due to the 
molecule itself. 



(a) 




(b) 




(c) 



pure molecule 

molecule with Au-cluster 




0.001 



E [eV] 

FIG. 3: (a) (Colour online) Structure for a molecule without 
cluster, (b) (Colour online) Sturcture for a molecule with a 
gold cluster on top. (c) Transmission function T(E) for both 
structures. Energy scale is relative to the HOMO-LUMO gap, 
such that E = corresponds to the middle of the gap. 

Therefore if transmission is already suppressed by the 
molecule (at all off-resonant energies), it can either be 
further reduced by off-resonant tunneling through the 
cluster, or it can (at best) be left unchanged by resonant 
tunneling through the cluster. Under no circumstances 
can transmission, once suppressed by the molecule, be 
afterwards increased by the cluster. This in turn means, 
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that metallic clusters can not give rise to additional peaks 
in the transmission spectrum. 



2. Mono-layer versus single molecule 



(a) 




i " ' 
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FIG. 4: (Colour online) The structure for one, two, three, 
and four molecules adsorbed within an Au-9 x 3 super-cell. 
This setup was used to test the sum rule. 

How do we expect the transmission function T l (E) 
for i periodically arranged molecules to look like? As 
long as the inter-molecular interactions are small (com- 
pared to the intra-molecular ones) the molecular levels 
of each molecule will not be significantly changed. Fur- 
thermore as the mono-layer consists of only one kind 
of molecule, all of them will have the same electronic 
structure. Therefore we expect each molecule to con- 
tribute the same amount to the transmission function: 
T n {E) := Y,i Tl ( E ) = nT 1 (E), where i runs over all n 
adsorbed molecules. 

We calculated the transmission function for n=l to 
4 molecules within an Au super-cell of size 9x3 (the 
structures are shown in Fig. |IJ. The distance between 
the molecules is chosen to be a multiple of the closest Au- 
Au separation a (d = 5.76A = 2a, with a = 2.88A). To 
our knowledge, the parameters of the PDFSAM mono- 
layer have never been determined experimentally, which 
is why we have to assume the above values. However 
STM studies^ and also theoretical calculations^ have 
been performed for alkanethiol mono-layers, and these 
parameters motivated our choice. 

Independent from the number of molecules present, the 
transmission functions have the same amount of peaks, 
at identical energetic positions (see Fig.[S]a). This result 
is also obtained for all larger distances of the molecules, 
where the inter- molecular interaction is even smaller. 
Furthermore, the sum rule is indeed fulfilled, as shown 
in Fig. 0(b). Each T l (E) is plotted against T 4 (E) for 
i G {1, 2, 3}. The calculated transmission values (300 dis- 
crete values each) clearly show a linear correlation. The 
straight lines are linear fits to the data, and their slope 
does very well agree with the theoretically expected value 
of o(n, m) = n/m. The deviation is below 6%, as can be 
seen in table [I] where we summarise all the fitted values 
for T n {E) = a{n, m) ■ T m (E). 

We conclude the following: a mono-layer, where the 
inter-molecular distance is large enough to not let inter- 
molecular interactions play a significant role, has the 
same number of distinct electronic levels as a single 
molecule. These levels are then highly degenerate. A 
CV-diagram will therefore have the same number of 
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a(3.4) = 0.7579 

a(2,4) = 0.51 17 

a(l,4) = 0.3649 




T*(E) 



FIG. 5: (a) Transmission function for one, two, three, and 
four PDI molecules, (b) When plotted against each other, the 
transmission functions reveal a linear relationship: T Z (E) = 
a{i,4)T i (E) (bottom). 



peaks. Only the net current will be increased by a factor 
a(n, m) compared to the single molecule case. 

The mere fact, that one deals with a mono-layer in- 
stead of a single molecule does not imply that the trans- 
mission function changes qualitatively. 



3. Influence of molecular clusters 

We have seen, that one does not observe additional 
peaks in the transmission function, as long as the inter- 
molecular influence is small. And this is the case for dis- 
tances which occur in typical SAM structures We 
now investigate cases, where the molecular interactions 
are not negligible. This occurs for example, when the pe- 
riodic structure of the mono-layer is perturbed by an ad- 
ditional molecule, such that a molecular cluster is formed. 
It is sufficient to study the transmission function of an 
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n 


m 


<x(n, m) 


-f- - 1 

n/m 


a 


1 


2 


0.5188 


3.76% 


0.014 


1 


3 


0.3500 


4.99% 


0.018 


1 


4 


0.2649 


5.96% 


0.021 


2 


3 


0.6753 


1.30% 


0.008 


2 


4 


0.5117 


2.34% 


0.015 


3 


4 


0.7579 


1.06% 


0.011 



TABLE I: The fitted values for o(n, m) together with their 
deviation from the theoretical value a(n,m) := n/m and 
a measure for the quality of the fit a, where a 2 := (N — 
I) -1 E {T n (Ei) - a{n,m)T m (Ei)f for N = 300 discrete en- 
ergy values T{Ei). 



isolated cluster only, because we have already seen that 
molecules in the periodic SAM arrangement do not influ- 
ence each other. The sum of the transmission function 
for the periodic SAM and the transmission function for 
the molecular cluster is, due to the sum rule, the total 
transmission function for defect and SAM. 

We study the influence of a shorter distance between 
two, three, and four molecules on the transmission spec- 
trum and relate it to the discrete energies of the iso- 
lated molecules. The molecules are now separated by 
d = 2.88A which corresponds to the Au-Au atom spac- 
ing. The atomic structure for this calculation is shown 
in Fig. |B](a), the resulting transmission functions in Fig. 
©(b) and (c). 

By reducing the molecular separation from d\ = 2 ■ 
2.88A = 5.76A to d 2 = 2.88A, the transmission func- 
tion qualitatively changes. The number of peaks roughly 
doubles. The new peak positions are different from the 
ones we have obtained in the previous calculations. And 
this time, the peak positions do depend on the number of 
molecules involved. This is an important point, because if 
there are several molecular clusters with different molec- 
ular distances, then they all give rise to peaks at different 
energy values. The resulting transmission function is the 
sum of the individual functions and will thus contain far 
more peaks, than the transmission function for the non 
perturbed periodic layer. 

We now show that the new peaks are a result of the 
increase in molecular interaction due to the decrease in 
spatial separation. For non-interacting molecules, the 
molecular energies are identical and therefore degenerate. 
An interaction between molecules breaks this degeneracy 
and therefore new energy levels occur. By performing a 
diagonalisation of the molecular Hamiltonian (without 
leads) one can determine the levels of the molecular clus- 
ter. 

In Fig. we have again plotted the transmission func- 
tion for three and four molecules, this time together with 
the discrete energy levels of the corresponding molecular 
cluster. The inset is identical to Fig. (c), while the 
plot itself is a magnification, to better resolve the dis- 
crete energy levels (which are shown as points along the 




-a o a 

E [eV] 




four molecules 

: I i i i I i i i I i i i I i l! 

-3 2 4 

E [eV] 

FIG. 6: (a) (Colour online) Two, three, and four molecules 
with a shorter inter-molecular distance, (b) The transmission 
functions for two and three molecules, (c) The transmission 
functions for three and four molecules. In contrast to all pre- 
vious cases, the peaks are shifted with respect to each other 
and there are also additional peaks. These changes are due 
to the increase in inter-molecular interaction, which alters the 
electronic levels. 

transmission function). Each of the transmission peaks 
is related to at least one discrete energy value. But not 
each energy value can be related to a peak in the trans- 
mission function. Why is that? The discrete energies can 
only give rise to new peaks in the transmission function, 
if they are not suppressed by a weak coupling to one of 
the leads. All levels which are not related to any peak 
belong to this category. If the position of the peak is 
shifted away from a corresponding energy level, then this 
is due to the coupling between molecules and leads. This 
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FIG. 7: Magnification of the transmission functions for three 
and four closely spaced molecules. Additionally the discrete 
energy levels of the system without leads are plotted as points 
along the transmission function. To each peak there belongs 
at least one discrete energy level. A detailed discussion is 
given in the main text. Inset: Transmission function (original 
scale) for three and four closely spaced molecules (identical 
to Fig. E| c). 

coupling is absent in the diagonalisation of the molec- 
ular Hamiltonian, but present in the calculation of the 
transmission function. 

Finally we show, that the additional peak structure in 
the transmission function for a scenario with an increased 
inter- molecular interaction gives rise to a number of steps 
in the TV-curve. Figure [3] contains an IV calculation 
for a molecular structure containing all three molecular 
clusters shown in Fig. El (a). In this calculation, the bias 
voltage Vb enters as a shift of the Fermi levels for source 
and drain lead: [i\ = ji% + sVb- The molecular energy 
has been set to E m = fix — 5E m — rjeV, where SE m is 
the zero bias displacement of the molecular levels and 
T] = 0.5, because of the symmetric coupling to the leads. 

Compared to the experiments*!^ the number of steps 
in the TV-curve is well reproduced by our calculation. 
The obtained current is at least one order of magnitude 
larger than the experimental values. 14 This is a phe- 
nomenon common to all theoretical methods based on 
the Landauer formulaAS^ A satisfactory explanation for 
this discrepancy as well as for the broad range of experi- 
mentally observed current values has not yet been found. 



the qualitative structure of the transmission function. 
Only if the distance between molecules becomes so small, 
that inter-molecular interactions are no longer negligible 
(which is below 6A in our case), does the transmission 
function undergo a qualitative change. Namely an addi- 
tional peak structure occurs. 

How does this finding compare to the experimental 
data? As we have pointed out in section [n] only in 
devices using molecules with two cyanide end-groups a 
more or less random peak structure was observed in the 
CV characteristic!*!* 1 ^. In other devices, molecules with 
at least one thiol end-group are typically used. These 
show significantly less peak structure. 

We therefore give the following interpretation: The 
thiol end-group is known to bind strongly to Au atoms. 
It is therefore likely, that thiol-based mono-layers sta- 
bly adsorb to gold leads. Resulting periodic structures 
are then robust against distortions. The conductance of 
such structures is proportional to the corresponding sin- 
gle molecule conductance, i.e. the number of molecules 
involved changes the absolute value of the current only, 
not the peak structure. 

The random like peak structure in devices, made up of 
cyanide based molecules suggests, that there are a some 
molecular clusters present in the mono-layer. These clus- 
ters might occur, because the binding of a cyanide end- 
group to Au is considerably weaker compared to that of 
a thiol end-group, and weaker binding results in a less 
robust periodic structure. 



V. DISCUSSION 



We have shown, that the peak structure of the trans- 
mission function is robust against changes in the number 
of adsorbed molecules, as long as the distance between 
molecules is considerably large id > 6A). And also does 
the exact shape of the top metallic lead not influence 
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APPENDIX A: CONNECTION BETWEEN 
EIGEN- AND CURRENT- VALUES 

The current properties of each channel can be related 
to the corresponding eigenvalue. We start from equation 



GD): 



WWjW) =(iP\W, 



=AIA 2 ^|W^|0>. 
6), i.e. Ai = A 2 , i.e. 



\WM) = 



Let's first consider |i/>) = 
lAKV'IWjl^). For each channel with eigenvalue |A| ^ 1 
one then must have (^W^l^) = 0, i.e. this channel does 
not itself carry any current. This is consistent with our 
terminology of an evanescent wave. If, however, |A,| = 1, 
then (ip\Wj is purely imaginary, because Wj is an anti- 
hcrmitian operator. We can therefore define the velocity 
of a propagating wave to be Vi := Im(?/>| Wj 

Now we consider the case of two different solutions 
^ |0) and define Ui )2 ■= (ip\Wj\(j)) . If their eigenval- 
ues do not satisfy A1A2 = 1, then the current between 
these two solutions is zero v\p = 0. So let's assume 
Ai = 1/A 2 . Because if |Ai| > 1 then | Aa j < 1, a current 
can flow between an evanescent left going wave and an 
evanescent right going wave. But if we restrict ourselves 
to solutions with finite amplitudes in a semi-infinite lead, 
then either the left or right going wave amplitude must 
be zero. Therefore evanescent waves do neither carry a 
current themselves nor do they exchange current with 
other channels, that is they do not at all contribute to 
the net current. 

Finally we are left with the case Ai = 1/A 2 *, with 
l^i I = IA2I = 1. This is equivalent to Ai = A2, i.e. 
the case of degenerate eigenvalues. Therefore propagat- 
ing waves to degenerate eigenvalues do exchange current. 
That in turn means, that the current of a superposition 
of two such waves does not necessarily equal the sum of 
the two individual currents, which is problematic as we 
want to express the total current as a sum of indepen- 
dent channels. However, the propagating and evanes- 
cent waves were obtained by diagonalising the propaga- 
tor P. This transformation is unique up to rotations 
in every degenerate eigenvalue subspace. Because W is 
anti-hermitian we can diagonalise these subspaces and 
the resulting diagonal elements will be purely imaginary. 
So the net current may be written as a summation over 
all the individual contributions of propagating channels, 
only if these subspace rotations are performed. 

Summarising we have shown that the transformation 
U diagonalising the propagator P (i.e U^ 1 PU) can be 
chosen such that the transformation U^WU of the cur- 
rent operator is diagonal in the subspace of propagating 
waves with purely imaginary diagonal elements. All the 
other diagonal entries are zero and the only non-zero non- 
diagonal elements belong to evanescent waves in opposite 
directions. 

APPENDIX B: CALCULATION OF THE 
SCATTERING MATRIX 



(Using this order for the coefficients it is straight for- 
ward to extend all formulas to the general case of more 
than two leads.) The indices 1 and 2 indicate source and 
drain lead surface layers, while the index is used for 
the molecular region. ti j2 are the coupling matrices from 
source/drain to the molecules. 

We now transform into the basis of incoming and out- 
going channels, i.e. we apply 



U = 



u 1 
u 2 

1 



with IP = 



from the right to equation (|B1(I : 



(H-ES)U 







A% A% t 
Bl B\ Bl B\ 



t 
1 



(B2) 



with 



At 



hiU- 



> - r M i Ui > A i > , 



and 



Bl = TjJ7>A>. 



The first and third column act on the surface layer of the 
incoming channels, the second and fourth act on outgo- 
ing ones, while the fifth column, acting on the molecular 
region, remains unchanged. 

The scattering matrix expresses the outgoing channel 
amplitudes in terms of the incoming ones. Therefore we 
split the matrix of equation l|B2[l into two parts, one con- 
taining the outgoing columns, the other one containing 
the incoming ones as well as the molecular column: 



Mn 



M in := 



A\ 4 

A\ 4 
B< Bl- M Q 



Ai 
A% 
Bl Bl 



and 



The first matrix Af out is a square matrix and by inverting 
it, we obtain the scattering matrix 



-M- u \ ■ M in . 



(B3) 



The part of the Hamiltonian containing the molecular 
region and its coupling to the leads can be written as 



(H-ES)\i>) = 



h x Mi r\ 
h 2 M 2 4 
n r 2 M 



|V) =0. (Bl) 
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